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ABSTRACT 

Classical design methods involved in magnetic bearings 
and magnetic suspension systems have always had their 
limitations. Because of this, the overall effectiveness of 
a design has always relied heavily on the skill and 
experience of the individual designer. This paper 
combines two approaches that have been developed to aid 
the accuracy and efficiency of magnetostatic design. The 
first approach integrates classical magnetic circuit theory 
with modem optimization theory to increase design 
efficiency. The second approach uses loss factors to 
increase the accuracy of classical magnetic circuit theory. 
As an example, an axial magnetic thrust bearing is 
designed for minimum power. 


NOMENCLATURE 1 

A g = Area of the pole pieces at the air gap 
B = Flux density 

B r = Remanence flux of neodymium iron boron 

F = Force 

J = Current density 

K a = Actuator loss factor 

K F = Flux leakage factor 

Kj= Coil mmf loss factor 

L g = Length of air gap 

L m = Permanent magnet thickness (on each pole piece) 
V c = Volume of the coil 
T| = Coil packing factor 
p 0 = Permeability of free space (4tix 1 O' 7 H/m) 
p = Resistivity of copper (2x1 O' 5 Q-m) 

0) = Frequency of occurrence of disturbance force 


INTRODUCTION 

One of the most important factors to consider in the 
design of an active magnetic bearing is its power 
efficiency. Reducing the power consumption will reduce 
not only the operational costs, but also the acquisition 
costs. If the required power can be reduced, smaller, and 
therefore cheaper amplifiers and power supplies can be 
used. Reduction of power also lessens winding cooling 
requirements. While applications such as space based 
systems and small pumps for medical applications have 
the most to gain by reducing power requirements, more 
typical commercial systems can also be significantly 
improved. 

Various approaches have been developed to aid the 
design of magnetic bearings. Parametric and systematic 
approaches have been developed to optimize the power to 
weight ratio of magnetic bearings [1] and to achieve 
minimum power [2]. These approaches are very useful in 
developing design strategies but do not typically yield 
specific solutions. 

A direct minimum power optimization method was 
developed in 1998 [3]. However, it optimized an ideal 
model with no flux leakage or mmf losses considered. It 
yielded helpful design trends but did not give design- 
specific solutions. This approach has now been extended 
with the incorporation of finite element analysis. 

Magnetic circuit theory is used to mathematically model 
the bearing. The weaknesses of classical circuit theory 
are overcome by including loss factors to account for flux 
leakage and magnetomotive force losses. The extended 
circuit theory is then incorporated into an optimization 
package. 


All equations are formulated in SI units 




The bearing geometry that is being optimized is shown in 
Fig. 1. The system is being optimized for minimum 
power. The bearing is designed to support a load of 2025 
N at an air gap of 4 mm. It is biased by a pair of 
neodymium iron boron permanent magnets placed on the 
poles of the stator. The remanent flux of the neodymium 
iron boron is 1.2 Tesla. The shaft clearance has a radius 
of 25 mm. The flux density in the iron is assumed to 
behave linearly up to a level of 1 Tesla. The maximum 
allowable current density of the windings is limited to 4 
A/mm 2 . The windings have a packing factor of 0.85. The 
overall volume of the coil is limited to 820 cm 3 . 



Figure 1: Axial Thrust Bearing Geometry 


MODELING 

Magnetic circuit theory was developed as an analogy to an 
electric circuit to model the behavior of a magnetic 
actuator. Several assumptions must be made for this 
analogy to hold: (1) there is no fringing or leakage flux, 
(2) there are no coil mmf losses, and (3) the permeability 
of the actuator core and return path is infinite. Applying 
magnetic circuit theory with these assumptions, the force 
applied to the rotating shaft by the magnetic actuator 
shown in Fig. 1 is [4], 


The losses due to the finite permeability of the stator and 
armature can be combined to form the actuator loss 
coefficient K a . The coil mmf loss can be accounted for by 
Kj. With the introduction of these terms into Eq. 1, the 
extended magnetic circuit model becomes, 


4 (K a L g + K F L m ) 2 


( 2 ) 


The three coefficients are unknown at this point but can 
be determined from finite element analysis. If the 
actuator has already been built, the values can be 
determined from experimental data. The values of K a and 
K F are determined by rearranging Eq. 2 and setting NI = 
0. 



Fpea is the force calculated by a finite element program. 
To accurately determine K a and K p , Fpp A must be 
determined at several different air gap distances. The 
more points that are used, the more accurate the solution. 
The validity of the solution is limited to the calculated 
range of air gap distances. The least squares solution for 
K a and K F can then be determined using the pseudo 
inverse. 



Once K a and K F have been determined, K, can be 
determined from Eq. 5. 


F = 





(i) 


The accuracy of Eq. 1 depends on how well the 
assumptions are met. If the assumptions are reasonably 
well satisfied, the accuracy of the model can be within a 
few percent. However, in common practice the error can 
be much larger. 


The accuracy can be greatly increased by relaxing these 
assumptions. This is done by introducing loss factors. 
Following the procedures outlined in [4], three factors can 
be introduced. The flux leakage around the magnetic 
circuit can be accounted for by the loss coefficient K F . 



2 B r L m 
0 NI 


The appropriate choice of NI used to solve for K, is 
important. If too small a value of NI is chosen, then the 
product KjNI ~0, so the value of K, becomes arbitrary. 
The most accurate method for determining K ( is to analyze 
several current levels near J max and solve for the least 
squares solution. However, choosing a single value of NI 
near the upper limit of the operating range will provide a 
suitable value of K,. 
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OPTIMIZATION 

It is desired to optimize an axial magnetic thrust bearing 
for minimum power. The bearing will be designed to 
provide a specified force at a given gap distance. The 
bearing will be able to handle disturbance forces without 
exceeding the maximum allowable current density of the 
windings, and without saturating the magnetic flux within 
the circuit. 

The power dissipated by a coil can be determined using 
the standard power equation, 

P = pr\J 2 V c (6) 


can be developed by substituting NI = T|A C J into Eq. 9 and 
solving for J. 


j Max, Min j 

r K a L g + K F L m > 

1 2B r L m 

1 ' 

v n4 j 

1 h„ti4 


The second inequality constraint requires that the flux 
density within the iron remains below the saturation level, 
1 Tesla, at every point within the circuit. To avoid 
contractions where the flux may saturate, the minimum 
cross sectional area of the iron was specified to equal the 
pole area at all points. 


At the operating point, the permanent magnets should 
provide a bias flux to carry the nominal load on the 
bearing. The power required by the coil will only be a 
factor when disturbance forces are present. For this 
particular case, it is assumed that the frequency of the 
disturbances follows a Gaussian distribution. 


cV 2 tc" 



The magnitude of the disturbance forces for this case have 
been set to ± 10% of the operating force. Once the 
optimization code settles on an actuator configuration, the 
current density required to provide the designated loads 
can be written as a function of force. With this, the 
objective function can be written as, 

F r 2 

f = pi\V c J 2 (F)e 2[ ° 1 dF 


0.2 F , 

*M V c-—[L J ( F i) e ‘ 


The performance characteristics of the actuator are 
controlled through the constraint equations. The bearing 
is being designed to support a specified load at a specified 
air gap. This introduces the equality constraint. 


H 0 A I KjNI + 


2 B r L m 
Ho 


4 (K a L g +K F L m ) 2 


There are three additional inequality constraints that must 
be added to the system. The first requires that the current 
density required to support the load at the high and low 
end of the design range remain below J max . The equation 


B Sat ~ 


Ho 


2 B r L m "j 

Ho J 


2 (K a L g +K F L n ) 


( 11 ) 


The final constraint requires the volume of the actuator to 
remain finite. For this case the maximum volume was 
limited to 820 cm 3 . 

Vmax ( 12 ) 

The coil geometry has been divided as shown in Fig. 2. 
Because the cross sectional area of the circuit was 
specified to equal the pole area, the outer wall thickness 
t, and the back wall thickness b, can be written in terms of 
the other variables. The design vector can therefore be 
reduced to 5 variables, 


x=\r w h L m J J (13) 



Figure 2: Description of Actuator Geometry 


The geometry of the actuator is not known at the start of 
the optimization process, so no loss factors can be 
calculated. An estimation of the final optimization space 
can be determined by setting all the loss coefficients to 1. 
Thus, the assumptions associated with classical circuit 
theory are active and it is found that the actual force 
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produced by the coil will typically be less than that 
predicted by the optimized solution. 

The ideal model is optimized and a solution is obtained. 
The actuator is then analyzed using finite element analysis 
and the loss factors are obtained. These loss factors are 
then introduced into the optimization code and the 
extended circuit theory model is optimized. Assuming 
that the addition of the loss factors does not significantly 
change the coil geometry, the optimum actuator 
configuration is then determined. If the geometry does 
change significantly, additional steps will be needed to 
determine precise values of K a ,K F , and Kj for the 
optimizer. 


IMPLEMENTATION OF EQUATIONS 

MATLAB’ s optimization toolbox version 5.3. was chosen 
to optimize the axial thrust bearing. The package uses the 
Sequential Quadratic Programming Method (SQP), which 
is ideally suited for nonlinear problems with nonlinear 
constraints. SQP uses the BFGS quasi-Newton method to 
approximate the Hessian of the Lagrangian function. The 
new Hessian is used to solve a Quadratic Programing 
subproblem in order to determine the best search direction 
and step length [5]. MATLAB also offers the advantage 
that the optimization toolbox is commercially available. 

The finite element modeling was performed using Vector 
Fields PC-OPERA 2D. Axial symmetry was used with the 
modified rxa potential [6]. The calculations were 
performed using quadratic elements with non-linear 
material properties. The iron was modeled using the 
standard material BH curve supplied by Vector Fields, 
and the Neodymium Iron Boron had a remanent flux of 
1.2 Tesla and a coercivity of 9.7xl0 5 A/m. 

There was some dependence upon the initial conditions to 
whether or not the optimizer would converge. When it 
did converge it was always to the same solution. When it 
did not converge, the initial conditions were randomly 
perturbed by +10% until the solution converged. To keep 
the units consistent, all lengths were implemented in the 
optimization code as meters. Accordingly, the current 
density was scaled to vary between 0 and 1 since the order 
of magnitude of r, w, h, and L m was 0.1 while the order of 
Jmax was 1 e6. 

Once the optimizer settled on a configuration six finite 
element runs were used to determine the loss factors. K a 
and K F were determined using five air gap lengths about 
the operating point, L g = 3.0, 3.5,4.0,4.5, and 5.0 mm. K ; 
was determined by setting NI = T)A c J max . 


DISCUSSION OF RESULTS 

The ideal model was optimized with the loss factors set to 
1 and a minimum power configuration was found. The 
resulting configuration and actual loss factors are given as 
Design-A in Table 1. The configuration delivered a zero 
power force of 1778 N, 12% below the desired value. 
The loss factors were then used in the optimization loop 
to develop the extended optimized model, Design-B. The 
addition of the constant loss factors only increased the 
zero power force to 1800 N, an improvement of less then 
2%. The ideal and extended models of Design-A, along 
with the extended model of Design-B, and the finite 
element results are shown in Fig. 3. 



Figure 3: Actuator Performance for Designs A and B 


The addition of the loss factors in the extended model of 
Design-A increased the accuracy of the formulation from 
15% to 0.3% when compared to the finite element results. 

Including the constant loss factors in the optimization 
loop had a negligible effect on the desired performance 
level. This is because the extended optimized design 
varied significantly from the ideal. The loss factors K a , 
K F , and K, for Design-B had an error of 10%, 2%, and 
1%, respectively, when compared to the actual values. 
This is due to the optimizer not being aware of the 
additional losses, due to the changes in the actuator 
geometry, during the optimization process. 

Design-B had a narrower pole face which increased flux 
leakage. It also had a thicker permanent magnet which 
increased the affective air gap reluctance. The increase in 
the conductor height, h, additionally increased the 
reluctance of the iron path. The full effect of these 
changes were not accounted for in the optimization 
process and resulted in a poor configuration. Thus, the 
addition of constant loss factors in the optimization loop 
does not provide a global minimum power configuration. 
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To dynamically include the loss factors in the 
optimization loop the loss factors must change with the 
geometry of the actuator and accurately represent the loss 
characteristics. Currently, a method is under development 
that will provide an estimation of the loss factors over a 
range of possible actuator geometries. This technique is 
still being developed and is the subject of future work. 
Using a preliminary approach to including estimated loss 
factors into the optimization loop resulted in Design-C, 
listed in Table 1. 


TABLE 1: Optimization Results and Loss Factors 



The performance results are shown in Figs. 4 and 5. The 
actuator geometry for Design-C is illustrated in Fig. 6. 



Air Gap (mm) 


Figure 4: Actuator Performance for Designs A, 
B, and C at NI = 0 


• Finite Element Results 



-1000 -500 0 500 1000 

NI (Amps) 

Figure 5: Actuator Performance for Designs A, 

B, and C at L g = 4mm. 



I* 26.6 mm 17.4 mW 

5.0 mm 

Figure 6: Optimized Actuator Configuration 


Comparing designs A and C reveals several trends. The 
pole area of Design-C increased which reduced the flux 
leakage. It also increased the volume of the permanent 
magnet which was necessary to compensate for the losses 
between the ideal model and the extended model. The 
permanent magnet volume increased while keeping the 
magnet thickness essentially the same, reducing the 
apparent gap reluctance. Also, the height of the conductor 
was reduced. This lessened the volume of the copper 
windings and reduced the reluctance of the iron path. 

Inclusion of the varying loss factors in the optimization 
loop resulted in a zero current force of 2046 N, 1% from 
desired. With the varying loss factors in the optimization 
loop, the actuator performance improved by 268 N at the 
zero current condition for the same size actuator. The 
accuracy of this model is evident by how well the Design- 
C model predicts the finite element results in Figs. 4 & 5. 

The optimization routine converged to a minimum power 
configuration. The solution can be checked by perturbing 
the design variables and analyzing the objective function. 
The solution is presumed to be the global minima within 
the accuracy of the estimated loss factors. For Design-C, 
the estimated loss factors K a , K F , and K ; varied from the 
exact values by 1.4%, 1.4%, and 3% respectively. 
Improvements to these values would slightly modify the 
optimum results but should not produce any significant 
changes. 

By reducing the reluctance losses and flux leakage at the 
gap, the actuator loss factor decreased by 3%. The flux 
leakage factor and coil mmf loss factor did not vary much 
between the two designs. This was because the slight 
variations between the two configurations did not 
dramatically change the overall actuator characteristics. 






Efficiency trends evident to the optimizer for the ideal 
case and the extended case are different. In the ideal 
model, the winding volume is minimized by increasing 
the height of the conductor and reducing the mean radius. 
This leads to a tall-narrow design that reduces power by 
reducing winding volume. In the extended model with the 
varying loss factors, the optimizer is aware that this design 
increases the flux loss and converges to a short-wide 
design. The same opposing efficiency trends are evident 
in the permanent magnet results. The tall-narrow design 
increased the permanent magnet volume by making the 
magnet thicker. In Design-C, the permanent magnet 
volume was increased by making the magnet wider. 


CONCLUSIONS 

This design method provides a unique approach for the 
design of magnetic bearings. Not only is the final design 
optimally suited to the specifications of the designer, but 
it also results in an accurate model if the finite element 
results are accurate. In the axial thrust bearing example 
covered in this paper, the extended circuit model had an 
accuracy of 99%. 

The method allows design constraints to be easily 
incorporated. Examples include constraints on flux 
density, current density, stiffness, volume, and geometric 
constraints. While in this case the bearing was optimized 
for minimum power, the objective statement could be 
easily modified to maximize or minimize force, stiffness, 
weight, volume, as well as many other possibilities. 

The use of constant loss factors to estimate losses in the 
optimization loop did not prove successful. The minimum 
power trends within the ideal and the extended models 
oppose each other. In order to optimize the realistic 
model, the variation of the loss factors with geometry 
must be accounted for. 

In the example, the combined use of the optimization code 
and variable loss factors increased the zero power force by 
268 N. The performance was improved with no increase 
in actuator volume and only a slight increase in permanent 
magnet volume. 
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